1 Introduction

The following text originally came from the Ghosh UZH UFO application

As the global climate is experiencing more heat, and less rainfall - it is reasonable to expect that the distributions of these variables (temperature and rainfall) are becoming more skewed and asymmetric towards the extreme values (see figure 1.1 below, which comes from the UZH UFO proposal, where it was figure 2). With the availability of more open access long-term databases, it is possible to address how different taxa respond at the community level.

Introduction figure

Figure 1.1: Introduction figure

As a preliminary work, I have already gathered long-term (median of 41 years) species-level abundance data for 2043 terrestrial and 716 aquatic communities. My recent result (manuscript in preparation) shows that the community stability is different for terrestrial and freshwater taxa and could be better explained considering the different strengths between pairwise species associations at the extremes, called community-level tail association, than with classic correlates of community stability studies (richness and variance ratio).

I will gather global data for annual temperature and rainfall from open access CHELSA database19 and ask how variability in temperature and precipitation would affect terrestrial taxa (birds, mammals, invertebrates, plants). For freshwater taxa (fish, phytoplankton, invertebrates), mostly temperature variability would be considered. In marine realm, sampling is spatially not consistent over the years20 and also very few long-term (>20yrs) data sampled compared to terrestrial and freshwater, thus I will focus only on the latter two realms. Also, species-level biomass (or body size) data will be gathered considering different generation times across taxa.

I will focus on community stability and will build a Bayesian model incorporating climatic factors (e.g., variability, skewness, range of maximum and minimum of temperature-distribution over the years etc.). While scientists studied thermophilization in the context of warming-related turnover in communities, no predictive model for community stability has been developed, to date, assessing the effect of extreme climatic events across taxa using a global database. This study will inform the current status of communities across multiple taxa facing climatic extremes and help prioritize conservation efforts (see Work Package 2).

I will gather annual climate data (mean, minimum, maximum for temperature, rainfall) and compute the variability, and the skewness of their distribution for the study period over which the community dynamics was studied. I will compute the richness (number of total species and dominant species that were present minimum 70% of the total years sampled), variance ratio, community level total tail association from pairwise synchrony as drivers. These drivers appeared as significant for explaining variation in community stability from my recent study (manuscript in preparation). I will compute the response variable community stability as the inverse of community-variability over the study period. Then, I will build a Bayesian model to see the effect of climate parameters, on the stability-driver relationships for different taxa.

2 Data

2.1 Data structure

Total 1947 community timeseries we have collected for the timespan 1979-2019. 4 taxa are considered - birds, fish, freshwater invertebrates, terrestrial invertebrates. Below is the summary of the datatable. Description of each column is given in README.txt

## 'data.frame':    1947 obs. of  55 variables:
##  $ source                  : chr  "BioTIME" "BioTIME" "BioTIME" "BioTIME" ...
##  $ STUDY_ID                : chr  "57" "229" "229" "229" ...
##  $ newsite                 : chr  "57" "STUDY_ID_229_LAT35.04016_LON-83.36127" "STUDY_ID_229_LAT35.11187_LON-83.39091" "STUDY_ID_229_LAT35.14137_LON-83.29577" ...
##  $ REALM                   : chr  "Freshwater" "Freshwater" "Freshwater" "Freshwater" ...
##  $ TAXA                    : chr  "fish" "fish" "fish" "fish" ...
##  $ ORGANISMS               : chr  "fish" "fish" "fish" "fish" ...
##  $ initR                   : int  76 24 30 32 25 37 22 30 26 31 ...
##  $ nsp                     : int  34 14 11 17 14 16 13 11 14 10 ...
##  $ nyr_used                : int  32 23 20 21 21 23 23 20 20 28 ...
##  $ startyr                 : int  1981 1990 1990 1991 1990 1990 1990 1995 1995 1979 ...
##  $ endyr                   : int  2012 2013 2013 2012 2013 2013 2013 2014 2014 2006 ...
##  $ nint                    : int  561 91 55 136 91 120 78 55 91 45 ...
##  $ nind                    : int  503 59 35 97 61 91 67 46 60 36 ...
##  $ npos                    : int  35 29 20 27 29 27 11 7 26 7 ...
##  $ nL                      : int  24 21 13 9 3 14 3 3 15 7 ...
##  $ nU                      : int  11 8 7 18 26 13 8 4 11 0 ...
##  $ nneg                    : int  23 3 0 12 1 2 0 2 5 2 ...
##  $ L                       : num  3.027 2.557 1.046 1.028 0.267 ...
##  $ U                       : num  -1.252 -0.633 -0.49 -2.547 -2.938 ...
##  $ f_nind                  : num  0.897 0.648 0.636 0.713 0.67 ...
##  $ f_nL                    : num  0.0428 0.2308 0.2364 0.0662 0.033 ...
##  $ f_nU                    : num  0.0196 0.0879 0.1273 0.1324 0.2857 ...
##  $ f_nneg                  : num  0.041 0.033 0 0.0882 0.011 ...
##  $ cvsq_real               : num  1.565 0.189 0.239 0.11 0.158 ...
##  $ cvsq_indep              : num  1.5001 0.0963 0.0879 0.0475 0.1099 ...
##  $ phi                     : num  1.04 1.96 2.72 2.33 1.44 ...
##  $ phi_LdM                 : num  0.454 0.552 0.388 0.442 0.542 ...
##  $ skw_real                : num  5.172 0.363 0.505 2.41 -0.296 ...
##  $ skw_indep               : num  5.064 0.612 0.737 0.841 -0.877 ...
##  $ phi_skw                 : num  1.021 0.593 0.685 2.866 0.337 ...
##  $ iCV                     : num  0.799 2.303 2.045 3.008 2.512 ...
##  $ iCValt                  : num  1.81 1.9 1.56 3.7 2.15 ...
##  $ LONGITUDE               : num  -89.5 -83.4 -83.4 -83.3 -83.5 ...
##  $ LATITUDE                : num  44 35 35.1 35.1 35.2 ...
##  $ t_med                   : num  2809 2865 2868 2864 2864 ...
##  $ t_skw                   : num  0.4702 0.1315 0.1711 -0.0306 0.0381 ...
##  $ t_var                   : num  0.00401 0.00228 0.00208 0.00233 0.00239 ...
##  $ t_med_celcius           : num  7.77 13.31 13.67 13.23 13.2 ...
##  $ t_skw_celcius           : num  0.4702 0.1315 0.1711 -0.0306 0.0381 ...
##  $ t_var_celcius           : num  0.145 0.0492 0.0438 0.0504 0.0518 ...
##  $ t.lm.slope              : num  0.341 0.258 0.106 0.365 0.2 ...
##  $ t.lm.slope.sig          : int  1 0 0 0 0 0 0 0 0 0 ...
##  $ t.sens.slope            : num  0.333 0.333 0.215 0.448 0.319 ...
##  $ t.sens.slope.sig        : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ t.lm.slope.celcius      : num  0.0341 0.0258 0.0106 0.0365 0.02 ...
##  $ t.lm.slope.sig.celcius  : int  1 0 0 0 0 0 0 0 0 0 ...
##  $ t.sens.slope.celcius    : num  0.0333 0.0333 0.0215 0.0448 0.0319 ...
##  $ t.sens.slope.sig.celcius: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ GiniSimpson             : num  0.817 0.57 0.92 0.819 0.521 ...
##  $ Simpson                 : num  0.142 0.152 0.555 0.256 0.138 ...
##  $ Shannon                 : num  0.67 0.512 0.863 0.695 0.471 ...
##  $ Heip                    : num  0.291 0.22 0.693 0.385 0.19 ...
##  $ McIntosh                : num  0.658 0.428 0.852 0.688 0.384 ...
##  $ SmithWilson             : num  0.415 0.35 0.627 0.315 0.325 ...
##  $ Pielou                  : num  0.19 0.194 0.36 0.245 0.178 ...

But see the below table which shows the sample size for each taxa and the datasource, we have very few sample size for terrestrial invertebrates. I feel it’s better to write a paper about north american birds vs european fish (atleast we have >500 datapoints for birds and fish). But I am open to other ideas. I don’t know which kind of data requirement we need for response diversity, but if we can have the body size or biomass (as trait) data then I can test the H0: whether response diversity increases the stability or influenced by temperature?

## # A tibble: 12 × 3
##    TAXA                      source              n
##    <chr>                     <chr>           <int>
##  1 birds                     BBS              1227
##  2 birds                     BioTIME            19
##  3 fish                      BioTIME            11
##  4 fish                      BioTIMEx           25
##  5 fish                      RivFishTIME       544
##  6 freshwater invertebrates  BioTIME            15
##  7 freshwater invertebrates  BioTIMEx            2
##  8 freshwater invertebrates  InsectRoel         78
##  9 freshwater invertebrates  SwissLakeZoo        5
## 10 freshwater invertebrates  Zooplankton2014     7
## 11 terrestrial invertebrates BioTIME            11
## 12 terrestrial invertebrates BioTIMEx            3

2.2 Sitemap for each taxa

3 Methods

I want to see how community stability-drivers relationship would affect by the changing environmental variable (annual temperature distribution). Temperature could vary in many ways (see 3.1). I am considering three aspects of environmental (temperature) timeseries here: median of annual temperatures (\(t_{med}\)) during the study periods, trend (\(t_{trend}\)) and skewness (\(t_{skw}\)) of annual temperature timeseries for a given community. My intuition is:

  • community stability would be affected by any of these temperature component either directly or indirectly via the drivers.

  • My hypothesis for direct effect is: warming environment (i.e., increasing any of these temperature components) should make the total productivity (community level abundance or biomass) more variable - so should decrease community stability.

  • For indirect effect via richness: increasing temperature (i.e. \(t_{med}\), \(t_{trend}\)) might increase or decrease species richness and diversity-stability relationship could be affected via portfolio effect, or via the factors which are affected by the richness. e.g. terrestrial plants, beetles, and vertebrates show declining richness with incresing temperature in the past studies. [***Shya: add citations for birds and fish particularly which shows increasing or decreasing richness due to temperature change].

  • For indirect effect via overall synchrony (variance ratio): Similarly, changing \(t_{med}\), \(t_{trend}\) should relate how species interactions get modified in a changing environment, i.e., affecting the synchrony level (variance ratio). For example, as temperature increases (or decreases) \(t_{med}\), community might loose some species, and might be dominated by fewer species which may or may not be synchronous, depending on the species-traits at that particular environment. If the community is exposing to warmer environment, then maybe warm adapted species will become dominant and species with similar traits can have higher synchrony. I am not sure, what will happen, but we can see how community-level response (avg of all species’response to warming for a given community) changes with temperature across the communities? [***Shya: add citations for birds and fish particularly which shows increasing or decreasing synchrony due to temperature change].[***Shya: we can test this with standardized correlation, shown later in Sec. 4.2.]

  • For indirect effect via tail-dependent synchrony (total tail asymmetry): Community level tail dependent synchrony (A) is the measure of how much species pair are synchronous in their extreme values (rare or common). It is reasonable to think that with increasing \(t_{trend}\), species would be more similar in their extreme values - continuous warming/cooling could make the species assembly similar - either species become simultaneously rare (if all could not adapt well in warmer condition) or abundant (due to release from competitive pressure from other species which could not adapt well in warmer environment but were dominant otherwise). On the other hand, a finite \(t_{skw}\) means extreme events (like heatwaves or cold winter) which should be related to species’ extremes, especially if species response is regulated by some environmental threshold-like phenomena. Now depending on the tolerance limit (or temperature threshold) of all the species in a given community, their synchronous response in the extremes could differ. For example, if all species have similar tolerance limit then an extreme heat or winter event would make them more synchronous in their tails - so skewness would increase tail dependent synchrony. But if the species in the community vary in their tolerance limits then an extreme event (skewed temperature timeseries) would decrease the tail dependent synchrony. [***Shya: maybe we can test this with STI-related concept]

Temperature timeseries figure

Figure 3.1: Temperature timeseries figure

Temperature timeseries figure with real data

Figure 3.2: Temperature timeseries figure with real data

3.1 Variables estimated and modelled

(Perhaps make this into a table.)

Let \(N_{i,t,s}\) be the abundance (sometimes it was biomass data when abundance data were not available) of species \(i\) at time \(t\) at site \(s\). Total abundance at time \(t\) at site \(s\) is \(N_{t,s} = \sum_{i=1}^{s} N_{t,s,i}\).

Community stability at site \(s\) was estimated as the inverse of the coefficient of temporal variation in total community abundance (when abundance info were not available, then biomass): \(TempStab_s = 1 / CV(N_{t,s}) = abs(median(N_{t,s})) / IQR(N_{t,s})\)

Species richness at site \(s\) was estimated as the number of total species (\(nsp\)) and dominant species that were present minimum 70% of the total years sampled (\(R\)).

Species evenness at site \(s\) was estimated as Smith-Wilson matrix.

Community variance ratio: a measure of synchrony, scaled between 0 to 1 (Loreau & Mazancourt).

Community level total tail association from pairwise synchrony: see BioDyn project, Figure 1.

Temperature median: Median of CHELSA-extracted annual temperature timeseries for the study years included in the analysis for each community.

Temperature trend: Monotonic trend of annual temperature timeseries (computed by non-parametric Sen’s method or parametric linear fit slope). I used the Sen’s slope in the path model, as non-parametric estimation has some advantage, see wikipedia, but it is very similar to linear slope (see 4.16).

Temperature skew: Skewness of CHELSA-extracted annual temperature timeseries for the study years included in the analysis for each community.

Temperature variability: Temperature variability for the community during the study period = IQR(annual temperature distribution for the study period)/abs(median(annual temperature)). [Note: in Celcius scale temperature t_med could be negative for some sites, so use absolute value always]

4 Results

4.1 Community stability exploration

Stability-diversity relationship for birds and fish.

Figure 4.1: Stability-diversity relationship for birds and fish.

4.1.1 Birds

Stability-diversity relationship for birds at different temperature levels

Figure 4.2: Stability-diversity relationship for birds at different temperature levels

Stability-temperature relationship for bird communities at different richness levels

Figure 4.3: Stability-temperature relationship for bird communities at different richness levels

Stability-synchrony relationship for birds at different temperature levels

Figure 4.4: Stability-synchrony relationship for birds at different temperature levels

Synchrony-temperature relationship (scatterplot)

Figure 4.5: Synchrony-temperature relationship (scatterplot)

Synchrony-temperature relationship (boxplot)

Figure 4.6: Synchrony-temperature relationship (boxplot)

Synchrony richness relationship.

Figure 4.7: Synchrony richness relationship.

Synchrony temperature relationship.

Figure 4.8: Synchrony temperature relationship.

4.1.1.1 Basic statistics for birds

Model of bird stability:

term estimate std.error statistic p.value
(Intercept) 3.6679 0.4566 8.0325 0.0000
nsp 0.0285 0.0108 2.6342 0.0085
t_med_celcius -0.1290 0.0375 -3.4425 0.0006
nsp:t_med_celcius 0.0027 0.0009 2.9832 0.0029

Model of bird synchrony:

term estimate std.error statistic p.value
(Intercept) 0.2954 0.0218 13.5237 0.0000
nsp -0.0027 0.0005 -5.2091 0.0000
t_med_celcius 0.0033 0.0018 1.8185 0.0692
nsp:t_med_celcius -0.0001 0.0000 -1.4679 0.1424

4.1.1.2 Bird conclusions

Bird communities display a positive richness stability relationship. This relationship is stronger at higher temperatures. Equally, high richness bird communities are more stable at higher temperatures, while low richness bird communities are less stable at higher temperatures.

There is some suggestion that this may be explained by synchrony, but the statistics show no strong associations of synchrony with \(t_{med}\).

4.1.2 Fish

Stability-diversity relationship at different temperature levels

Figure 4.9: Stability-diversity relationship at different temperature levels

Stability-temperature relationship at different richness levels

Figure 4.10: Stability-temperature relationship at different richness levels

Stability-synchrony relationship at different temperature levels

Figure 4.11: Stability-synchrony relationship at different temperature levels

Synchrony-temperature relationship (scatterplot)

Figure 4.12: Synchrony-temperature relationship (scatterplot)

Synchrony-temperature relationship (boxplot)

Figure 4.13: Synchrony-temperature relationship (boxplot)

Synchrony richness relationship.

Figure 4.14: Synchrony richness relationship.

Synchrony temperature relationship.

Figure 4.15: Synchrony temperature relationship.

4.1.2.1 Fish basic statistics

Model of fish stability:

term estimate std.error statistic p.value
(Intercept) -0.3087 0.1559 -1.9800 0.0482
log2(nsp) 0.1873 0.0886 2.1133 0.0350
t_med_celcius 0.0638 0.0186 3.4343 0.0006
log2(nsp):t_med_celcius -0.0223 0.0087 -2.5625 0.0106

Model of fish synchrony:

term estimate std.error statistic p.value
(Intercept) -0.1910 0.1212 -1.5752 0.1158
log2(nsp) -0.3236 0.0689 -4.6968 0.0000
t_med_celcius -0.0293 0.0144 -2.0267 0.0432
log2(nsp):t_med_celcius 0.0075 0.0068 1.1036 0.2702

4.1.2.2 Fish conclusions

Fish communities display a positive richness stability relationship at low temperature, and a negative one at higher temperatures. This is the opposite of the interaction pattern for birds, where the relationship became more positive when temperature was higher. Equally, high richness fish communities are slightly less stable at higher temperatures, while low richness fish communities are more stable at higher temperatures (again the opposite to the bird patterns).

Not sure, at present, how much can be explained by synchrony, but the statistics show no strong associations of synchrony with t_med or richness.

Which pattern / finding are we trying to explain with the following two graphs?

Distribution of temperature trend estimated by non-parametric Sen's slope, and parametric linear fit slope. Colored points are significant Sen's slope (green: birds, blue: fish).

Figure 4.16: Distribution of temperature trend estimated by non-parametric Sen’s slope, and parametric linear fit slope. Colored points are significant Sen’s slope (green: birds, blue: fish).

## 
##  Freshwater Terrestrial 
##   0.3344828   0.2471910
Histogram plot for trends, both taxa.

Figure 4.17: Histogram plot for trends, both taxa.

4.2 Explanations

So, from the exploratory plots we can see: at higher temperature positive stability-diversity relationship becomes stronger for birds but for fish it becomes weaker. Also fish becomes more asynchronous with increasing temperature. So, why does that happen? to find this we could explore how much the bird species and fish species are consistent to temperature change across all communities.

The cue is: if fish species are not much consistent in their response to warming and vary across sites, that means you cannot make a conclusion that they would become similar with changing temperature. On another note, bird species should be more consistent towards warming if their is no change in their synchrony level across communities. Another possibility could be with changing temperature you might loose some species (its not just number of individuals, it will selectively prefer few species with better fitness), and then the communities will be dominated by few species with similar traits (so increasing synchrony). we will test this below.

From the above plots, we can see birds are showing consistent response-distribution across all temperature change, i.e., in either end of temperature spectrum (low or high end). That’s why the synchrony level remains similar for birds. But for fish, warming increases the richness (addition of new species), and as fish species now become more variable in response to temperature sensitivity (trait-variation), they show more asynchrony compared to low temperature scenario where only few species exists (see smaller circle size on the map for lowT,<50%CI) and show similar traits (so more synchrony). Note: when I show this to Frank, he commented on how much robust is the pattern for fish at low T as there are only few species existed across 145 sites - so it also depends on how we considered the lowT-highT communities. I set beyond 50% CI of temperature range as low/high. Even if I decrease that to 30% CI, still very few species found in low T sites (15 sp across 203 sites: 80% >0, 20% <0 line).

To further explore this idea: we collected traits data for birds and fish species used in the analysis. For fish-traits, I will use body length measurements, for bird-traits I will use HWI (Hand-wing index). From below figures: at high T, birds have slightly less dispersal ability (lower HWI), but richness is more or less uniformly spread at either temperature range. For fish, at lowT, few large species exists with similar traits (remember the previous histogram plot 90-10) showing higher synchrony, as temperature increases addition of new small fishes in the community (maybe better environment for them to exist in that temperature rather than too cold water) makes them asynchronous with more trait variation (histogram plot 66-34).

When I showed this to Blake, he was not convinced by the idea to split the data into two: low/high based on t_med (to him this temperature difference is more on latitudinal differences as shown in the map), and same species can exist in both communities - so why changing t_med should change the synchrony level for fish? and getting different bodysize fish from low/high t_med (fewer big fish in lowT and many smaller fish in highT) is not explaining why big fish should be more synchronous - is it because of fewer species (richness) or because bigger fish abundance change needs more time - not on annual scale?

So, I thought to make a plot of how community-level average response traits (average of standardised correlation between species abundance with t_med timeseries across sites) changes with increasing temperature (t_med)? For fish, it should decrease with increasing t_med, whereas for birds it should be a flat relationship.

Possible explanation:

Response variation with temperature

Figure 4.18: Response variation with temperature

Now, we will do a path analysis for a simplistic mixed effect model to see the environmental effects on community stability for both taxa.

## 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |==================                                                    |  25%
  |                                                                            
  |===================================                                   |  50%
  |                                                                            
  |====================================================                  |  75%
  |                                                                            
  |======================================================================| 100%
## 
## Structural Equation Model of psem_birds 
## 
## Call:
##   VR ~ R + E + MedianT + TrendT + SkewT
##   R ~ MedianT + TrendT + VarT
##   E ~ R + MedianT + TrendT + VarT
##   A ~ R + E + VR + MedianT + TrendT + SkewT
##   stability ~ R + E + VR + A + (R:MedianT) + MedianT + SkewT + TrendT + VarT
## 
##     AIC
##  14407.878
## 
## ---
## Tests of directed separation:
## 
##    Independ.Claim Test.Type       DF Crit.Value P.Value 
##   R ~ SkewT + ...      coef 1238.736     0.0072  0.9326 
##   E ~ SkewT + ...      coef 1239.471     3.1590  0.0758 
##   VR ~ VarT + ...      coef 1012.496     0.3846  0.5353 
##    A ~ VarT + ...      coef 1179.623     0.5959  0.4403 
## 
## --
## Global goodness-of-fit:
## 
## Chi-Squared = NA with P-value = NA and on 4 degrees of freedom
## Fisher's C = 8.191 with P-value = 0.415 and on 8 degrees of freedom
## 
## ---
## Coefficients:
## 
##    Response Predictor Estimate Std.Error        DF Crit.Value P.Value
##          VR         R  -0.3093    0.0322  840.4521    91.2304  0.0000
##          VR         E  -0.2142    0.0300 1061.4845    50.6442  0.0000
##          VR   MedianT   0.0219    0.0443  194.1824     0.2406  0.6243
##          VR    TrendT  -0.0294    0.0379  585.8074     0.5949  0.4408
##          VR     SkewT  -0.0165    0.0312  948.0008     0.2783  0.5979
##           R   MedianT  -0.1054    0.0532  649.3517     3.8844  0.0492
##           R    TrendT  -0.0004    0.0342 1231.4095     0.0002  0.9896
##           R      VarT  -0.1142    0.0267  822.8100    18.2258  0.0000
##           E         R  -0.0770    0.0315 1233.4093     5.9583  0.0148
##           E   MedianT   0.1246    0.0567  489.4629     4.7691  0.0294
##           E    TrendT   0.0846    0.0378 1188.9556     4.9779  0.0259
##           E      VarT   0.0189    0.0290  794.5467     0.4209  0.5167
##           A         R   0.7475    0.0224  464.0606  1099.9131  0.0000
##           A         E   0.1588    0.0213  716.9356    54.8898  0.0000
##           A        VR   0.6043    0.0209 1188.2246   833.3884  0.0000
##           A   MedianT  -0.0315    0.0244  140.0482     1.6356  0.2030
##           A    TrendT   0.0102    0.0237  260.7704     0.1807  0.6711
##           A     SkewT  -0.0543    0.0209  447.0396     6.6445  0.0103
##   stability         R   0.1463    0.0370  910.8877    15.4818  0.0001
##   stability         E  -0.0437    0.0260  897.3854     2.7891  0.0953
##   stability        VR  -0.5549    0.0318 1176.7268   303.5390  0.0000
##   stability         A  -0.1533    0.0330 1233.6583    21.4960  0.0000
##   stability   MedianT  -0.0310    0.0378  164.4429     0.6625  0.4169
##   stability     SkewT   0.0116    0.0254  583.8267     0.2063  0.6499
##   stability    TrendT   0.0346    0.0305  401.5987     1.2718  0.2601
##   stability      VarT  -0.0188    0.0241 1154.9180     0.6087  0.4354
##   stability R:MedianT   0.0478    0.0258  515.5481     3.3748  0.0668
##   Std.Estimate    
##        -0.3093 ***
##        -0.2142 ***
##         0.0219    
##        -0.0294    
##        -0.0165    
##        -0.1054   *
##        -0.0004    
##        -0.1142 ***
##        -0.0770   *
##         0.1246   *
##         0.0846   *
##         0.0189    
##         0.7475 ***
##         0.1588 ***
##         0.6043 ***
##        -0.0315    
##         0.0102    
##        -0.0543   *
##         0.1463 ***
##        -0.0437    
##        -0.5549 ***
##        -0.1533 ***
##        -0.0310    
##         0.0116    
##         0.0346    
##        -0.0188    
##         0.0544    
## 
##   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05
## 
## ---
## Individual R-squared:
## 
##    Response method Marginal Conditional
##          VR   none     0.15        0.32
##           R   none     0.01        0.65
##           E   none     0.03        0.57
##           A   none     0.58        0.60
##   stability   none     0.42        0.47
## 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |==================                                                    |  25%
  |                                                                            
  |===================================                                   |  50%
  |                                                                            
  |====================================================                  |  75%
  |                                                                            
  |======================================================================| 100%
## 
## Structural Equation Model of psem_fish 
## 
## Call:
##   VR ~ R + E + MedianT + TrendT + SkewT
##   R ~ MedianT + TrendT + VarT
##   E ~ R + MedianT + TrendT + VarT
##   A ~ R + E + VR + MedianT + TrendT + SkewT
##   stability ~ R + E + VR + A + (R:MedianT) + MedianT + SkewT + TrendT + VarT
## 
##     AIC
##  6310.489
## 
## ---
## Tests of directed separation:
## 
##    Independ.Claim Test.Type       DF Crit.Value P.Value 
##   R ~ SkewT + ...      coef 414.2328     3.7572  0.0533 
##   E ~ SkewT + ...      coef 236.9543     1.0553  0.3053 
##   VR ~ VarT + ...      coef 572.6102     0.7037  0.4019 
##    A ~ VarT + ...      coef 303.1369     0.3107  0.5777 
## 
## --
## Global goodness-of-fit:
## 
## Chi-Squared = NA with P-value = NA and on 4 degrees of freedom
## Fisher's C = 11.158 with P-value = 0.193 and on 8 degrees of freedom
## 
## ---
## Coefficients:
## 
##    Response Predictor Estimate Std.Error       DF Crit.Value P.Value
##          VR         R  -0.5503    0.0601 376.3677    81.8501  0.0000
##          VR         E  -0.2924    0.0401 567.0933    52.6589  0.0000
##          VR   MedianT  -0.2475    0.0588 116.3912    17.3370  0.0001
##          VR    TrendT   0.0674    0.0395 438.6255     2.8705  0.0909
##          VR     SkewT  -0.0939    0.0496 193.0384     3.5161  0.0623
##           R   MedianT   0.1327    0.0521 208.9389     6.4142  0.0121
##           R    TrendT  -0.0099    0.0284 575.8529     0.1214  0.7276
##           R      VarT  -0.0045    0.0260 515.8896     0.0296  0.8636
##           E         R  -0.4363    0.0607 452.1873    50.7083  0.0000
##           E   MedianT  -0.1376    0.0628 143.6693     4.7304  0.0313
##           E    TrendT   0.0440    0.0417 497.3481     1.1011  0.2945
##           E      VarT  -0.0272    0.0407 568.8198     0.4419  0.5065
##           A         R   1.0133    0.0323  59.1753   889.8806  0.0000
##           A         E   0.1100    0.0236 383.3190    21.1445  0.0000
##           A        VR   0.1107    0.0242 514.2745    20.6184  0.0000
##           A   MedianT  -0.1029    0.0270  54.6353    13.9146  0.0005
##           A    TrendT   0.0247    0.0213 251.1697     1.3024  0.2549
##           A     SkewT  -0.0046    0.0243 101.7171     0.0350  0.8521
##   stability         R   0.1803    0.1212 550.4539     2.1973  0.1388
##   stability         E  -0.0105    0.0481 568.6917     0.0474  0.8277
##   stability        VR  -0.3071    0.0473 562.8586    41.9415  0.0000
##   stability         A   0.0801    0.0777 546.5000     1.0605  0.3036
##   stability   MedianT   0.0745    0.0788 160.2071     0.8839  0.3485
##   stability     SkewT  -0.1454    0.0614 287.8530     5.5276  0.0194
##   stability    TrendT  -0.0001    0.0459 533.9515     0.0000  0.9983
##   stability      VarT   0.0379    0.0439 544.6907     0.7450  0.3884
##   stability R:MedianT  -0.2427    0.0708 568.6605    11.6996  0.0007
##   Std.Estimate    
##        -0.5503 ***
##        -0.2924 ***
##        -0.2475 ***
##         0.0674    
##        -0.0939    
##         0.1327   *
##        -0.0099    
##        -0.0045    
##        -0.4363 ***
##        -0.1376   *
##         0.0440    
##        -0.0272    
##         1.0133 ***
##         0.1100 ***
##         0.1107 ***
##        -0.1029 ***
##         0.0247    
##        -0.0046    
##         0.1803    
##        -0.0105    
##        -0.3071 ***
##         0.0801    
##         0.0745    
##        -0.1454   *
##        -0.0001    
##         0.0379    
##        -0.2089 ***
## 
##   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05
## 
## ---
## Individual R-squared:
## 
##    Response method Marginal Conditional
##          VR   none     0.39        0.49
##           R   none     0.04        0.53
##           E   none     0.25        0.42
##           A   none     0.78        0.78
##   stability   none     0.14        0.40

5 Discussion